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We propose a dynamical system that captures changes to the network centrality of 
nodes as external interest in those nodes vary. We derive this system by adding time- 
dependent teleportation to the PageRank score. The result is not a single set of impor- 
tance scores, but rather a time-dependent set. These can be converted into ranked lists 
in a variety of ways, for instance, by taking the largest change in the importance score. 
For an interesting class of the dynamic teleportation functions, we derive closed form so- 
lutions for the dynamic PageRank vector. The magnitude of the deviation from a static 
PageRank vector is given by a PageRank problem with complex-valued teleportation 
parameters. Moreover, these dynamical systems are easy to evaluate. We demonstrate 
the utility of dynamic teleportation on both the article graph of Wikipedia, where the 
external interest information is given by the number of hourly visitors to each page, and 
the Twitter social network, where external interest is the number of tweets per month. 
For these problems, we show that using information from the dynamical system helps 
improve a prediction task and identify trends in the data. 

1 Introduction 

The PageRank vector of a directed graph is the stationary distribution of a Markovian 
random surfer. At a node, the random surfer either 

1. transitions to a new node from the set of out-edges, or 

2. does something else (e.g. leaves the graph and then randomly returns) [Page et al., 
1998, Langville and Meyer, 2006]. 

The probability that the surfer performs the first action is known as the damping param- 
eter in PageRank, denoted a. The second action is called teleporting and is modeled by 
the surfer picking a node at random according to a distribution called the teleportation 
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distribution vector or personalization vector. This PageRank Markov chain always has 
a unique stationary distribution for any < a < 1. In this paper, we focus on the tele- 
portation distribution vector v and study how changing teleportation behavior manifests 
itself in a dynamical system formulation of PageRank. 

To proceed further, we need to formalize the PageRank model. Let A be the adjacency 
matrix for a graph where Aij denotes an edge from node i to node j. To avoid a 
proliferation of transposes, we define P as the transposed transition matrix for a random- 
walk on a graph: 

Pj i = probability of transitioning from node i to node j. 

Hence, the matrix P is column- stochastic instead of row-stochastic, which is the standard 
in probability theory. Throughout this manuscript, we utilize uniform random-walks on 
a graph, in which case P = where £> is a diagonal matrix with the degree of 

each node on the diagonal. However, none of the theory is restricted to this type of 
random walk and any column-stochastic matrix will do. If any nodes have no out-links, 
we assume that they are adjusted in one of the standard ways [Boldi et al., 2007]. Let 
V be a teleportation distribution vector such that Wi > and '^^Vi = 1. This vector 
models where the surfer will transition when "doing something else." The PageRank 
Markov chain then has the transition matrix: 

aP -I- (1 - Q;)ve^. 

While finding the stationary distribution of a Markov chain usually involves computing 
an eigenvector or solving a singular linear system, the PageRank chain has a particularly 
simple form for the stationary distribution vector x: 

{I-aP)yL= (l-a)v. 

Sensitivity of PageRank with respect to v is fairly well understood. Langville and 
Meyer [2006] devote a section to determining the Jacobian of the PageRank vector with 
respect to v. The choice of v is often best guided by an application specific measure. By 
setting V = e^, that is, the zth canonical basis vector: 



■ith row, 



PageRank computes a highly localized diffusion that is known to produce empirically 
meaningful clusters and theoretically supported clusters [Andersen et al., 2006, Tong 
et al., 2006]. By choosing v based on a set of known-to-be- interesting nodes, PageRank 
will compute an expanded set of interesting nodes [Gyongyi et al., 2004, Singh et al.. 
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Static PageRank Dynamic PageRank 

Figure 1 - At left, we have PageRank with static teleportation. At each step, the teleportation is 
to each node with uniform probabihty 1/5. At right, we have the PageRank model with dynamic 
teleportation. In this case, the teleportation distribution (illustrated in red) changes with time. 
Thus, the upper nodes are teleported to more during the middle time regime. In both cases, 
the graph is fixed. We study the effect of such dynamic teleportation on the PageRank scores. 

2007]. Yet, in all of these cases, v is chosen once for the graph application or particular 
problem. 

In the original motivation of PageRank [Page et al., 1998], the distribution v should 
model how users behave on the web when they don't click a link. When this intuition is 
applied to a site like Wikipedia, this suggests that the teleportation function should vary 
as particular topics become interesting. For instance, in our experiments (section 5), 
we examine the number of page views for each Wikipedia article during a period where 
a major earthquake occurred. Suddenly, page views to earthquake spike - presumably 
as people are searching for that phrase. We wish to include this behavior into our 
PageRank model to understand what is now important in light of a radically different 
behavior. One option would be to recompute a new PageRank vector given the observed 
teleporting behavior at the current time. Our proposal for a dynamical system is another 
alternative. That is, we define a new model where teleportation is the time-dependent 
function: 

vit). 

At each time t, v(t) is a probability distribution of where the random walk teleports. 
Figure 1 illustrates this model. We return to a comparison between this approach and 
solving PageRank systems in section 4. 

The dynamical system we propose is a generalization of PageRank in the sense that 
if v(t) is a constant function in time, then we converge to the standard PageRank vector 
(theorem 4). Additionally, we can analyze the dynamical PageRank function for some 
simple oscillatory teleportation functions v(t). Bounding the deviation of these oscil- 
latory PageRank values from the static PageRank vector involves solving a PageRank 
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problem with complex teleportation [Horn and Scrra-Capizzano, 2007, Constantino and 
Glcich, 2010]. This result is, perhaps, the first non-analytical use of PageRank with a 
complex teleportation parameter. 

In our new dynamical system, we do not compute a single ranking vector as others 
have done with time-dependent rankings [Grindrod et al., 2011], rather we compute a 
time-dependent ranking function x(t), the dynamic PageRank vector at time t, from 
which we can extract different static rankings (section 2.4). There are two complications 
with using empirically measured data. First, we must choose a time-scale for our ODE 
based on the period of our page- view data (section 2.5). Put a bit informally, we must 
pick the time-unit for our ODE - it is not dimensionless. We analytically show that 
some choices of the time-scale amount to solving the PageRank system for each change 
in the teleportation vector. Second, we also investigate smoothing the measured page 
view data (section 2.6). To compute this dependent ranking function x(<) we discuss 
ordinary differential equation (ODE) integrators in section 3. 

We discuss the impact of these choices on two problems: page views from Wikipedia 
and a retweet network from Twitter. We also investigate how the rankings extracted 
from our methods differ from those extracted by other static ranking measurements. We 
can use these rankings for a few interesting applications. Adding the dynamic PageRank 
scores to a prediction task decreases the average error (section 6.1) for Twitter. Clustering 
the dynamic PageRank scores yields many of the standard time-series features in social 
networks (section 6.2). Finally, using Granger causality testing on the dynamic PageRank 
scores helps us find a set of interesting links in the graph (section 6.3). 

We make our code and data available in the spirit of reproducible research: 

http; / /www. cs .purdue . edu/homes/dgleich/ codes/dynsyspr-im 



2 PageRank with time-dependent teleportation 

We begin our discussion by summarizing the notation introduced thus far in Table 1 . 

In order to incorporate changes in the teleportation into a new model for PageRank, 
we begin by reformulating the standard PageRank algorithm in terms of changes to the 
PageRank values for each page. This step allows us to state PageRank as a dynamical 
system, in which case we can easily incorporate changes into the vector. 

The standard PageRank algorithm is the power method for the PageRank Markov 
chain [Langville and Meyer, 2006]. After simplifying this iteration by assuming that 
e^x = 1, it becomes: 

^ aPx^^) + (1 - a)v. 

In fact, this iteration is equivalent to the Richardson iteration for the PageRank linear 
system (J — aP)x = (1 — a)v. This fact is relevant because the Richardson iteration is 
usually defined: 



x(fe+i) = x(^) +0J (1 - a)v - (/ - aP)x 



(fe) 
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Table 1 - Summary of notation. Matrices are bold, upright roman letters; vectors are bold, 
lowercase roman letters; and scalars are unbolded roman or greek letters. Indexed elements are 
vectors if they are bolded, or scalars if unbolded. 



I 


the imaginary number 


n 


number of nodes in a graph 


e 


the vector of all ones 


P 


column stochastic matrix 


a 


damping parameter in PageRank 


V 


teleportation distribution vector 


X 


solution to the PageRank computation; (7 — aP)'K = (1 — q)v 


x(t) 


solution to the dynamic PageRank computation for time t 


v(t) 


a teleportation distribution vector at time t 


c 


the cumulative rank function 


r 


the variance rank function 


d 


the difference rank function 


Vfc 


the teleportation distribution for the fcth observed page-views vector 


e 


decay parameter for time-series smoothing 


s 


the time-scale of the dynamical system 


^max 


the last time of the dynamical system 



For = 1, we have: 

AxC^) = xC^+i) - x^^) = aPx^^) + (1 - a)v - x^^) = (1 - a)v - (7 - aP)x('^). 

Thus, changes in the PageRank values at a node evolve based on the value (1 — a)v — 
(7 — q;P)x*^'^^ . We reinterpret this update as a continuous time dynamical system: 

x'(<) = (1 -a)v- (7-aP)x(t). (1) 

To define the PageRank problem with time-dependent teleportation, we make v{t) a 
function of time. 

DEFINITION 1 The dynamic PageRank model with time- dependent teleportation is the 
solution of 

x'(t) = (l-a)v(t)-(7-aP)x(t) (2) 

where x(0) is a probability distribution vector and v{t) is a probability distribution vector 
for all t. 

In the dynamic PageRank model, the PageRank values x(t) may not "settle" or 
converge to some fixed vector x. We see this as a feature of the new model as we plan to 
utilize information from the evolution and changes in the PageRank values. For instance, 
in section 2.4, we discuss various functions of x(i) that define a rank. Next, we state the 
solution of the problem. 
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LEMMA 2 The solution of the dynamical system: 



:K.'{t) = (1 - a)v(i) - (/ - aP)x{t) 

is ^ 

x(t) = exp[-(J - aP)t]x(0) + (1 - a) / exp[-(J - aP){t - T)]v(r) dr. 

Jo 

This result is found in standard texts on dynamical systems, for example [Berman et al., 
1989]. 

Given this solution, let us quickly verify a few properties of this system: 

LEMMA 3 The solution of a dynamical PagcRank system x(i) is a probability distribution 
(x(t) > and e^x(t) ^ 1) for all t. 

Proof The model requires that x(0) is a probability distribution. Thus, x(0) > and 
e^x(O) — 1. Assuming that the sum of x(t) is 1, then the sum of the derivative x'(f) 
is as a quick calculation shows. The closed form solution above is also nonnegative 
because the matrix exp[— (/ — aP)] = exp(Q;P) exp(— 1) > and both x(0) and v{t) are 
non-negative for all t. (This property is known as exponential non-negativity and it is 
another property of Af-matrices such as / — aP [Herman et al., 1989].) ■ 

2.1 A generalization of PageRank 

This closed form solution can be used to solve a version of the dynamic problem that 
reduces to the PageRank problem with static teleportation. If v{t) = v is constant with 
respect to time, then 

/ exp[-(/ - aP){t - t)]v(t) dr = {I - aPy^v - exp[-(/ - aP)t]{I - aPy^v. 
Jo 

Hence, for constant v(t): 

x(<) = cxp[-(/ - aP)i](x(0) - x) + X, 

where x is the solution to static PageRank: (/ — aP)x = (1 — q;)v. Because all the 
eigenvalues of —(7 — aP) are less than 0, the matrix exponential terms disappear in a 
sufficiently long time horizon. Thus, when v(i) = v, nothing has changed. We recover 
the original PageRank vector x as the steady-state solution: 

lim x(<) = X the PageRank vector. 

t—^oo 

This derivation shows that dynamic teleportation PageRank is a generalization of the 
PageRank vector. We summarize this discussion as: 
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THEOREM 4 PageRank with time- dependent teleportation is a generalization of PageR- 
ank. If w{t) — w , then the solution of the ordinary differential equation: 

x'(<) = (1 - a)v - (/ - aP)x{t) 

converges to the PageRank vector 

[I - aP)yi = (1 - a)^r 

as < — >■ oo. 

2.2 Choosing the initial condition 

There are three natural choices for the initial condition x(0). The first choice is the 
uniform vector x(0) — ^e. The second choice is the initial teleportation vector x(0) = 
v(0). And the third choice is the solution of the PageRank problem for the initial 
teleportation vector (/ — aP)x(O) = (1 — a)v(O). We recommend either of the latter two 
choices in order to generalize the properties of PageRank. Note that if x(0) is chosen 
to solve the PageRank system for v(0), then x(t) ~ x for all t is the solution of the 
PageRank dynamical system with constant teleportation (theorem 4). 

2.3 PageRank with fluctuating interest 

One of the advantages of the PageRank dynamical system is that we can study problems 
analytically. We now do so with the following teleportation function, or forcing function 
as it would be called in the dynamical systems literature: 

At) = \Y.^,{'^o,{t + {j-iy-^) + i), 

where Vj is a teleportation vector. Here, the idea is that Vj represents the propensity 
of people to visit certain nodes at different times. To be concrete, we might have vi 
correspond to news websites that are visited more frequently during the morning, V2 
correspond to websites visited at work, and V3 correspond to websites visited during 
the evening. This function has all the required properties that we need to be a valid 
teleportation function. With the risk of being overly formal, we'll state these as a lemma. 

LEMMA 5 Let k > 2. Let vi,...,Vfc be probability distribution vectors. The time- 
dependent teleportation function 

^(.t) = lj2^,[cosit + ij-l)^) + l), 

satisfies the both properties: 
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1. v{t) > for all t, and 
2- ElLi = 1 for all t 

Proof The first property follows directly because the minimum value of the cosine func- 
tion is —1, and thus, w{t) is always non-negative. The second property is also straight- 
forward. Note that 

n k k 

J2v{t)i = 1 + J2^os{t + (j - l)f ) = 1 + 5]Re{exp(zt + t{j - l)f )}• 

i=l j=l j=l 

Let rj{t) = exp{tt + i{j — 1)^). For t = 0, these terms express the k roots of unity. For 
any other t, we simply rotate these roots. Thus we have rj (t) — for any t because 
the sum of the roots of unity is if fc > 2. The second property now follows because the 
sum of the real component is still zero. ■ 

For this function, we can solve for the steady-state solution analytically. 

LEMMA 6 Let k > 2, < a < 1, P be column-stochastic, vi,...,Vfc be probability 
distribution vectors, and 

^(^) = ^ E^^ H^^ + (j' - + l) = lvcos{t + f) + ive. 



where V ■ 
of 



Vl,...,Vfe 



3 = 1 

2tt 



and fj = (j — 1)^ j = 1, . . . ,k. Then the steady state solution 
x'(t) = (1 - a)v{t) -{I- aP)x{t) 



IS 

x{t) = X -|- Re{sexp(jt)} 
where x is the solution of the static PageRank problem 

(I-aP)x = (1 - a)lve 
k 

and s is the solution of the static PageRank problem with complex teleportation 

(7-^P)s = (l-a)^Fcxp(*f). 

Proof This proof is mostly a derivation of the expression for the solution by guessing the 
form. First note that if 

x(t)=x + y(t) 

then 

y'{t) = (1 - a)^Vcos{t + f) - (/ - aP)y{t). 



8 



That is, we've removed the constant term from the teleportation function by looking 
at solutions centered around the static PageRank solution. To find the steady-state 
solution, we look at the complex-phasor problem: 

•L'it) = (1 - a)iyexp(«i + if) - (I-aP)'z{t) 

where y(t) = Re{z(t)}. Suppose that z(i) = sexp(it). Then: 

z'(t) = zsexp(ii) — (1 — a)i exp(zt) exp(if) — (/ — Q!P)s exp(it). 

The statement of s in the theorem is exactly the solution after canceling the phasor 
exp(it). We now have to show that this solution is well-defined. PageRank with a complex 
teleportation parameter 7 exists for any column-stochastic P if I7I < 1 (see Horn and 
Serra-Capizzano [2007], Constantine and Gleich [2010]). For the problem defining s, 
— a/{l + i) and I7I = aj \f2. Thus, such a vector s always exists. ■ 

We conclude with an example of this theorem. Consider a four node graph with 
adjacency matrix and transition matrix: 



A = 



0010 
0010 
0101 
1100 



and 



00 0.5 

0.5 0.5 
110 
0.5 



Let Vj = ej for j = 1, ... ,4. That is, interest oscillates between all four nodes in the 
graph in a regular fashion. We show the evolution of the dynamical system for 20 time- 
units in Figure 2. This evolution quickly converges to the oscillators predicted by the 
lemma. In the interest of simplifying the plot, we do not show the exact curves as they 
are visually indistinguishable from those plotted for i > 4. By solving the complex valued 
PageRank to compute s, we can compute the magnitude of the fluctuation: 



s| = 0.0216 0.0261 0.0122 0.0235 



This vector accurately captures the magnitude of these fluctuations. 




■ Page 1 

■ Page 2 

■ Page 3 

■ Page 4 



Figure 2 - The dashed lines 
represent the average PageR- 
ank vector computed for the 
teleportation vectors. The 
curves show the evolution of 
the PageRank dynamic sys- 
tem for this example of tele- 
portation. We see that the 
dynamic PageRanks fluctuate 
about their average PageRank 
vectors. Lemma 6 predicts the 
magnitude of the fluctuation. 
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2.4 Ranking from Time-Series 



The above equations provide a time-series of dynamic PageRank vectors for the nodes, 
denoted formaUy as x(t), < t < imax- AppHcations, however, often want a single score, 
or small set of scores, to characterize sets of interesting nodes. There are a few ways in 
which these time series give rise to scores. Many of these methods were explained by 
O'Madadhain and Smyth [2005] in the context of ranking sequences of vectors. Having a 
variety of different scores derived from the same data frequently helps when using these 
scores as features in a prediction or learning task [Bccchetti et al., 2008, Constantine and 
Gleich, 2010]. 

Transient Rank. We call the instantaneous values of x(t) a node's transient rank. 
This score gives the importance of a node at a particular time. 

Summary, Variance, & Cumulative Rank. Any summary function s of the time 
series, such as the integral, average, minimum, maximum, or variance, is a single score 
that encompasses the entire interval [0,tniax]- We utilize the cumulative rank, c and 
variance rank, r in the forthcoming experiments: 




Difference Rank. A node's difference rank is the difference between its maximum and 
minimum rank over all time, or a limited time window: 

d = max[x(t)l — minfxft)] dw = maxx(f) — minx(i). 

t t ^ ^ " tew tew 

Nodes with high difference rank should reflect important events that occurred within the 
range [0,tinax] or the time window W. We suggest using a window W that omits the 
initial convergence region of the evolution. In the context of Figure 2, we'd set W to 
be [4,20] to approximate the vector |s| numerically. In Section 6 and figure 7, we see 
examples of how current news stories arise as articles with high difference rank. 



2.5 Modeling activity 

In the next two sections of our introduction to the dynamic teleportation PageRank 
model, we discuss how to incorporate empirically measured activity into the model. Let 
Pi, . . . ,Pf. he k observed vectors of activity for a website. In the cases we examine below, 
these activity vectors measure page views per hour on Wikipedia and the number of 
tweets per month on Twitter. We normalize each of them into teleportation distributions, 
and conceptually think of the collection of vectors as a matrix 



Vl, . 



V= Vl,, 



, Vfc 
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Let e{i) be a functional form representing the vector e^. The time-dependent teleporta- 
tion vector we create from this data is: 

v{t) = Fe(fioor{t} + 1) = Vfioo,{(}+i. 

For this choice, the time-units of our dynamical system are given by the time-unit of the 
original measurements. Other choices are possible too. Consider: 

Wsit) ^ ye(floor{t/s} + 1) = Vfloor{t/4+i- 

If s > 1, then time in the dynamical system slows down. If s < 1, then time accelerates. 
Thus, we call s the time-scale of the system. Note that 

x(sj), j = 0, 1,... 

represents the same effective time-point for any time-scale. Thus, when we wish to 
compare different time-scales s, we examine the solution at such scaled points. 

In the experimental evaluation, the parameter s plays an important role. We illustrate 
its effect in Figure 3(a) for a small subnetwork extracted from Wikipedia. As we discuss 
further in section 3, for large values of s, then v(t) looks constant for long periods of time, 
and hence x(t) begins to converge to the PageRank vector for the current, and effectively 
static, teleportation vector. Thus, we also plot the converged PageRank vectors as a step 
function. We see that as s increases, the lines converge to these step functions, but for 
s = 1 and s = 2, they behave differently. 

2.6 Smoothing empirical activity 

So far, we defined a time-dependent that v(t) changes at fixed intervals based on empiri- 
cally measured data. A better idea is to smooth out these "jumps" using an exponentially 
weighted moving average. As a continuous time function, this yields: 

v'{t;9) = 9v{t) - 0v{t;0). 

To understand why this smooths the sequence, consider an implicit Euler approximation: 

1 hf) 

This update can be written more simply as: 

v(<;^)- 7v(t) +a_7)v(i_M), 

new data old data 

where 7 = yq^. When v{t) changes at fixed intervals, then v(t; 0) will slowly change. 
If 9 is small then v{t; 9) changes slowly. We recover the "jump" changes in v(t) in the 
limit — >■ 00. 

The effect of 9 is shown in Figure 3(b). Note that we quickly recover behavior that is 
effectively the same as using jumps in v(<) {9 = 1, 10). So we only expect changes with 
smoothing for < 1. 
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(a) timescale s 



(b) smoothing 8 



— a = 0.5---a = 0.85 - - a = 0.99 




(c) damping parameter a 

Figure 3 - The evolution of PageRank values for one node due to the dynamical teleportation. 
The horizontal axis is time [0, 20], and the vertical axis runs between [0.01,0.014]. In figure (a), 
a = 0.85, and we vary the time-scale parameter (section 2.5) with no smoothing. The solid 
dark line corresponds to the step function of solving PageRank exactly at each change in the 
teleportation vector. All samples are taken from the same effective time-points as discussed in 
the section. In figure (b), we vary the smoothing (section 2.6) of the teleportation vectors with 
s = 2, and a — 0.85. In figure (c), we vary a with s — 2 and no smoothing. We used the ode45 
function in Matlab, a Runge-Kutta method, to evolve the system. 

2.7 Choosing the teleportation factor 

Picking a even for static PageRank problems is challenging, see Gleich et al. [2010] and Con- 
stantine and Gleich [2010] for some discussion. In this manuscript, we do not perform 
any systematic study of the effects of a beyond Figure 3(c). This simple experiment 
shows one surprising feature. Common wisdom for choosing a in the static case suggests 
that as a approaches 1, the vector becomes more sensitive. For the dynamic telepor- 
tation setting, however, the opposite is true. Small values of a produce solutions that 
more closely reflect the teleportation vector ~ the quantity that is changing - whereas 
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large values of a reflect the graph structure, which is invariant with time. Hence, with 
dynamic teleportation, using a small value of a is the sensitive setting. Note that this 
observation is a straightforward conclusion from the equations of the dynamic vector: 

x'(i) = (1 - a)v(i) + aPTi{t) - x(t) 

so a small implies a larger change due to v(t). Nevertheless, we found it surprising in 
light of the existing literature. 

3 Methods for dynamic PageRank 

In order to compute the time-sequence of PageRank values x(t), we can evolve the 
dynamical system (2) using any standard method - usually called an integrator. We 
discuss both the forward Euler method and a Runge-Kutta method next. Both methods, 
and indeed, the vast majority of dynamical system integrators only require a means to 
evaluate the derivative of the system at a time t given x(f). For PageRank with dynamic 
teleportation, this corresponds to computing: 

x'(i) = (1 - a)v(i) - (I - aP)x{t). 

The dominant cost in evaluating x'(i) is the matrix vector product Px. For the explicit 
methods we explore, all of the other work is linear in the number of nodes, and hence, 
these methods easily scale to large networks. Both of these methods may also be used 
in a distributed setting if a distributed matrix-vector product is available. 

3.1 Forward Euler 

We first discuss the forward Euler method. This method lacks high accuracy, but is fast 
and straightforward. Forward Euler approximates the derivative with a first order Taylor 
approximation: 

_ x(t + h)- x(0 
X Itj ~ 

and then uses that approximation to estimate the value at a short time-step in the future: 

x{t + h)= x(t) + h[{l- a)w{t) - (I- aP)x{t)] . 

This update is the original Richardson iteration with h — lo. We present the forward 
Euler method as a formal algorithm in Figure 4 in order to highlight a comparison with 
the power and Richardson method. That is, the forward Euler method is simply running 
a power method, but changing the vector v at every iteration. However, we derived this 
method based on evolving (2). Thus, by studying the relationship between (2) and the 
algorithm in Figure 4, we can understand the underlying problem solved by changing the 
teleportation vector while running the power method. 
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Input: 

a graph G — (V, E) and a procedure to compute Px for this graph 
a maximum time imax 

a function to return v(i) for any <t < tmax 
a damping parameter a 
a time-step h 

Output: X where the kth column of X is x(0 + kh) for aU 1 < fc < tmax/^ (or any 
desired subset of these values) 

x(0) -s— v(0) (or any other desired initial condition) 
while t < tinax — h do 

x{t + h) ^ x(t) + h [(1 - a)v{t) - (/ - aP)x{t)] 

X{:,k) ^x{t + h) 

t^t + h: k^k + 1 

Figure 4 - The forward Euler method for evolving the dynamical system: x'(t) = (1 — a)v{t) — 
(7 — aP)x(t). The resulting procedure looks remarkably similar to the standard Richardson 
iteration to compute a PageRank vector. One key difference is that there is no notion of 
convergence. 

Long time-scales. Using the forward Euler method, we can analyze the situation 
with a large time-scale parameter s. Consider an arbitrary x(0), a — 0.85, s ~ 100, 
h — 1, and no smoothing. In this case, then the forward Euler method will run the 
Richardson iteration for 100 times before observing the change in v(i) at t — 100. The 
difference between x(/c) and the exact PageRank solution for this temporarily static v(t) 
is ||x(fc) — x|ji < 2a''. For k > 50, this difference is small. Thus, a large s and no 
smoothing corresponds to solving the PageRank problem for each change in v. 

Stability. The forward Euler method with timestep h is stable if the eigenvalues of the 
matrix —h{I — aP) are within distance 1 of the point —1. The eigenvalues of P are all 
between —1 and 1 because it is a stochastic matrix, and so this is stable for any h < 

3.2 Runge-Kutta 

Runge-Kutta [Rungc, 1895, Kutta, 1901] numerical schemes are some of the most well- 
known and most used. They achieve far greater accuracy than the simple forward Euler 
method, at the expense of a greater number of evaluations of the function x'(i) at each 
step. We use the implementations of Runge-Kutta methods available in the Matlab ODE 
suite [Shampine and Reiclielt, 1997]. The step-size is adapted automatically based on a 
local error estimate, and the solution can be evaluated at any desired point in time. The 
stability region for Runge-Kutta includes the region for forward Euler, so these methods 
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are stable. These methods are also fast. To integrate the system for Wikipedia with over 
4 million vertices and 60 million edges, it took between 300-600 seconds, depending on 
the parameters. 

3.3 Maintaining interpretability 

Based on the theory of the dynamic teleportation system, we expect that x{t) > and 
e-^x(i) = 1 for all time. Although this property should be true of the computed solution, 
we often find that the sum diverges from one. Consequently, for our experiments, we 
include a correcting term: 

x'(t) = (1 - a)v(i) - (77 - aP)x(t) 

where 7 = (1 — a)e'^v{t) + ae^x(t). Note that 7 = 1 if the x{t) has sum exactly 1. If 
e^x(t) is slightly different from one, then the correction with 7 ensures that e-^x'(t) — 
numerically. Similar issues arise in computing static PageRank [Wills and Ipsen, 2009], 
although the additional computation in the Runge-Kutta methods exacerbates the prob- 
lem. 

4 Related work 

Note that we previously studied this idea in a conference paper [Rossi and Gleich, 2012]. 
These ideas have been significantly refined for this manuscript. 

The relationship between dynamical systems and classical iterative methods has been 
utilized by Embree and Lehoucq [2009] to study eigenvalue solvers. It was also noted 
in an early paper by Tsaparas [2004] that there is a relationship between the PageRank 
and HITS algorithms and dynamical systems. 

In the past, others studied PageRank approximations on graph streams [Das Sarma 
et al., 2008]. More recently, Bahmani et al. [2012] studied how accurately an evolving 
PageRank method could estimate the true PageRank of an evolving graph that is accessed 
only via a crawler. The method used here solved each PageRank problem exactly for 
the current estimate of the underlying graph. A detailed study of how PageRank values 
evolve during a web-crawl was done by Boldi et al. [2005]. In the future, we plan to study 
dynamic graphs via similar ideas. 

As explained in section 3 and figure 4, our proposed method is related to chang- 
ing the teleportation vector in the power method as its being computed. Bianchini et 
al. [Bianchini et al., 2005] noted that the power method would still converge if either the 
graph or the vector v changed during the method, albeit to a new solution given by the 
new vector or graph. Our method capitalizes on a closely related idea and we utilize the 
intermediate quantities explicitly. Another related idea is the Online Page Importance 
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Computation (OPIC) [Abitcboul et al., 2003], which integrates a PageRank-Uke compu- 
tation with a crawUng process. The method does nothing special if a node has changed 
when it is crawled again. 

While we described PageRank in terms of a random-surfer model, another character- 
ization of PageRank is that it is a sum of damped transitions: 

oo 

x = (l-a)^(aP)'=v. 

These transitions are a type of probabilistic walk and Grindrod et al. [Grindrod et al., 
2011] introduced the related notion of dynamic walks for dynamic graphs. We can 
interpret these dynamic walks as a backward Euler approximation to the dynamical 
system: 

x'(t) = aA{t)yi{t) x(0) = e 

with time-step h = l and A is a time-dependent adjacency matrix. This relationship sug- 
gests that there may be a range of interesting models between our dynamic teleportation 
model and existing evolving graph models. 

Outside of the context of web-ranking, O'Madadhain and Smyth propose Even- 
tRank [O'Madadhain and Smyth, 2005], a method of ranking nodes in dynamic graphs, 
that uses the PageRank propagation equations for a sequence of graphs. We utilize the 
same idea but place it within the context of a continuous dynamical system. In the con- 
text of popularity dynamics [Ratkicwicz et al., 2010], our method captures how changes 
in external interest influence the popularity of nodes and the nodes linked to these nodes 
in an implicit fashion. Our work is also related to modeling human dynamics, namely, 
how humans change their behavior when exposed to rapidly changing or unfamiliar con- 
ditions [Bagrow et al., 2011]. In one instance, our method shows the important topics 
and ideas relevant to humans before and after one of the largest Australian Earthquakes 
(figure 7). 

In closing, we wish to note that our proposed method does not involve updating the 
PageRank vector, a related problem which has received considerable attention [Chien 
et al., 2004, Langville and Meyer, 2004]. Nor is it related to tensor methods for dynamic 
graph data [Sun et al., 2006, Dunlavy et al., 2011]. 

5 Examples of dynamic teleportation 

We now use dynamic teleportation to investigate page view patterns on Wikipedia and 
user activity on Twitter. In the following experiments, unless otherwise noted, we set 
s = 1, a = 0.85, do not use smoothing ("0 = c»"), and use the ode45 method from 
Matlab to evolve the system. We study this model on two datasets. 
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5.1 Datasets 



We provide some basic statistics of the Wikipedia and Twitter datasets in Table 2. For 
Wikipedia, the time unit for s = 1 is an hour, and for Twitter, it is one month. 

Wikipedia Article Graph and Hourly page views. Wikipedia provides access to 
copies of its database [Various, 2009]. We downloaded a copy of its database on March 
6th, 2009 and extracted an article-by-article link graph, where an article is a page in the 
main Wikipedia namespace, a category page, or a portal page. All other pages and links 
were removed. See Gleich et al. [2007] for more information. 

Wikipedia also provides hourly page views for each page [Various, 2011]. These are 
the number of times a page was viewed for a given hour. These are not unique visits. 
We downloaded the raw page counts and matched the corresponding page counts to the 
pages in the Wikipedia graph. We used the page counts starting from March 6, 2009 
and moving forward in time. Although it would seem like measuring page views would 
correspond to measuring x(t) instead of v(t), one of our earlier studies showed that users 
hardly ever follow links on Wikipedia [Gleich et al., 2010]. Thus, wc can interpret these 
page views as a reasonable measure of external interest in Wikipedia pages. 

Twitter Social Network and Monthly Tweet Rates. We use a follower graph 
generated by starting with a few seed users and crawling follows links from 2008. We 
extract the user tweets over time from 2008—2009. A tweet is represented as a tuple (user, 
time, tweet). Using the set of tweets, we construct a sequence of vectors to represents 
the number of tweets for a given month. 



Table 2 - Dataset Properties. The page views or tweets is denoted as p. 



Dataset 


Nodes 


Edges 


^max 


Period 


Average pi 


Max Pi 


WIKIPEDIA 


4,143,840 


72,718,664 


48 


hours 


1.4243 


353,799 


TWITTER 


465,022 


835,424 


6 


months 


0.5569 


1056 



5.2 Rankings from transient scores 

First, we evaluate the rankings from dynamic PageRank using the intersection similarity 
measure [Boldi, 2005]. Given two vectors x and y, the intersection similarity metric at 
k is the average symmetric difference over the top-j sets for each j < k. li Xk and yt 
are the top-A: sets for x and y, then 



isimfc(x,y) = -}_^^— 
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(c) Average page views 
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Figure 5 - Intersection similarity of rankings derived from dynamic PageRank. We compute the 
intersection similarity of the difference, variance, and cumulative rankings given by dynamic 
PageRank and compare these with the rankings given by the in-degree, average page views, 
static PageRank with uniform teleportation, and static PageRank with average page views as 
the teleportation vector. For dynamic PageRank, we set the initial value x(0) to be the solution 
of the static PageRank system which uses v(0) as the teleportation vector. 



where A is the symmetric set- difference operation. Identical vectors have an intersection 
similarity of 0. 

For the Wikipedia graph, Figure 5 shows the similarity profile comparing a few rank- 
ing measures from dynamic PageRank to reasonable baselines. In particular, we compare 
d, r, c (from §2.4) to indegree, average page views, static PageRank with uniform telepor- 
tation, and static PageRank using average page views as the teleportation vector. The 
results suggest that dynamic PageRank is different from the other measures, even for 
small values of k. In particular, combining the external influence with the graph appears 
to produce something new. The only exception is in Fig. 5(d) where the cumulative rank 
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is shown to give a similar ordering to static PageRank using average page views as the 
teleportation. 

5.3 Difference ranks 

Figure 6 and figure 7 show the time-series of the top 100 pages by the difl^erence measure 
for Wikipedia with s = 1 and s ~ 4 without smoothing. Many of these pages reveal the 
ability of dynamic PageRank to mesh the network structure with changes in external 
interest. For instance, in figure 7, we find pages related to an Australian earthquake 
(43, 84, 82), the "recently" released movie "Watchmen" (98, 23-24), a famous musician 
that died (2, 75), recent "American Idol" gossip (34, 63), a remembrance of Eve Carson 
from a contestant on "American Idol" (88, 96, 34), news about the murder of a Harry 
Potter actor (60), and the Skittles social media mishap (94). These results demonstrate 
the effectiveness of the dynamic PageRank to identify interesting pages that pertain to 
external interest. The influence of the graph results in the promotion of pages such as 
Richter magnitude (84). That page was not in the top 200 from page views. 
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Figure 6 - The top- 100 Wikipedia pages that fluctuate the most as determined by the difi'erence ranking from our dynamic PageRank 
approach. The blue curves are the transient scores, and the red-curves are the raw page view data. The horizontal axis is 24 hours, and the 
vertical axes are normalized to show the range of the data. 
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Figure 7 - We again plot the top 100 pages from the difTerence rank, now using s = 4. The blue curves are the transient scores, and the 
red-curves are the raw page view data. The horizontal a^xis is 24 hours, and the vertical axes are normalized to show the range of the data. 
This choice gives a similar ordering as our previous forward Euler iteration method from Rossi and Gleich [2012]. Note the large change 
between the transient scores for "MainPage" between this figure and figure 6. 



6 Applications of time-dependent telportation 



This section explores the opportunity of using Dynamic PageRank for a variety of apph- 
cations outside of the context of ranking. 

6.1 Predicting future page views & tweets 

We begin by studying how well the dynamical system can predict the future. Formally, 
given a lagged time-series Pt_u,, P(_i, [Ahmed et al., 2010], the goal is to predict 
the future value p^+i (actual page views or number of tweets). This type of temporal 
prediction task has many applications, such as actively adapting caches in large database 
systems, or dynamically recommending pages. 

We performed one-step ahead predictions {t + 1) using linear regression. That is, we 
learn a model of the form; 

f (t - 1) i{t-2) ... f (t - w) 1 b w p(t) 

where w is the window-size, and f (•) is either page views or both page views and transient 
scores. After fitting b, the model predicts p(t + 1) as 

f{t) f(i-l) ••• f{t-w + l)] h 

We use the symmetric Mean Absolute Percentage Error (sMAPE) [Ahmed et al., 2010] 
measure to evaluate the prediction: 

sMAPE= ' ^ 1^*-^*" 



\T\^^ipt+Pt)/2- 

This relative error measure averages all the relative prediction errors over all the time- 
steps. We then average it over nodes. 

We study two predictive modes. The base model uses only the time-series of page 
views or tweets to predict the future page views or tweets. The dynamic teleportation 
model uses both the transient scores with smooting and page views to predict the future 
page views (or tweets). 

We evaluate these models for prediction on stationary and non- stationary time-series. 
Informally, a time-series is weakly stationary if it has properties (mean and covari- 
ance) similar to that of the time-shifted time-series. We consider the top and bottom 
10,000 nodes from the difference ranking as nodes that are approximately non-stationary 
(volatile) and stationary (stable), respectively. Table 3 compares the predictions of the 
models across time for non-stationary and stationary prediction tasks. Our findings 
indicate that the Dynamic PageRank time-series provides valuable information for fore- 
casting future tweet rates; however, it adds little (if any) accuracy in forecasting future 
page views on Wikipedia. 
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Table 3 - The ratio between the base model and the model with dynamic teleportation scores 
with s = 1,2,6, and oo, for three smoothing parameters. (Here, s — oo corresponds to solving 
the PageRank problem exactly for each change in teleportation.) If this ratio is less than 1, 
then the model with the dynamic teleportation scores improves the prediction performance. We 
also distinguish between prediction problems with highly volatile nodes (non-stationary) and 
nodes with relatively stable behavior (stationary) . The results show a much stronger benefit for 
Twitter than for Wikipedia 

Dataset Type 8 Error Ratio 

s (timcscalc) 







1 


2 


6 


oo 


TWITTER stationary 


0.01 


0.635 


0.929 


0.913 


0.996 




0.50 


0.636 


0.735 


0.854 


0.939 




1.00 


0.522 


0.562 


0.710 


0.963 


non-stationary 


0.01 


0.461 


0.841 


1.001 


0.992 




0.50 


0.261 


0.608 


0.585 


0.929 




1.00 


0.137 


0.605 


0.617 


0.918 


WIKIPEDIA stationary 


0.01 


0.978 


0.991 


0.989 


0.978 




0.50 


1.140 


1.130 


1.004 


0.990 




1.00 


1.084 


0.976 


1.010 


0.990 


non-stationary 


0.01 


0.968 


1.011 


0.968 


1.004 




0.50 


1.218 


0.994 


1.030 


1.031 




1.00 


1.241 


0.996 


0.957 


0.998 



For Twitter, the dynamic teleportation model improves predictions the most with the 
non-stationary nodes. The diffusion of activity captured by the model allows our model 
to detect, early on, when the external interest of vertices will change, before that change 
becomes apparent in the external interest of the vertices. This is easiest to detect when 
there is a large sudden change in external interest of a neighboring vertex. 

6.2 Clustering transient score trends 

Identifying vertices with similar time-series is important for modeling and understanding 
large collections of multivariate time-series. We now group vertices according to their 
transient scores. By using the difference rank measure d for s = 4, we cluster the top 
5,000 vertices using k-means with fc = 5, repeat the clustering 2,000 times, and take the 
minimum distance clustering identified. 

The cluster centroids are temporal patterns, and the main patterns in the dynamic 
PageRanks are visualized in figure 8(a). Pattern 2 represents European-centric behavior, 
whereas the others correspond to spikes or unusual events occurring within the dynamic 
PageRank system. Figure 8(b) plots the 20 closest vertices matching the patterns above. 
A few pages from the five groups are consistent with our previously discussed results 
from figure 7. One such unusual event is related to the death of a famous musician/actor 
from the Philippines (see pages 1-20). The pages from the third cluster (41-60) are 
related to "American Idol" and other TV shows/actors. Also some of the pages from 
the fourth cluster relate to Bernard Madoff (63, 66, 67, 70, 73), six days before he plead 
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Figure 8 - Vertices with similar dynamical properties are grouped together. The visualization re- 
veals the important dynamic patterns (spikes, trends) present from March 6th, 2009 in our large 
collection of time-series from Wikipedia. For each hour, we sample twice from the continuous 
function x(t) and utilize these intermediate values in the clustering. 



guilty in the largest financial fraud in U.S. history. This grouping reveals many of the 
standard patterns in time-series such as spikes and increasing/decreasing trends [Yang 
and Leskovcc, 2011]. 
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6.3 Towards causal link relationships 



In this section, we use Granger causality tests [Granger, 1969] on the collection of tran- 
sient scores to attempt to understand which links are most important. The Granger 
causality model, briefly described below, ought to identify a causal relationship between 
the time-series of any two vertices connected by a directed edge. This is because there 
is a causal relationship between their time-series in our dynamical system. However, 
due to the impact of the time-dependent teleportation, only some of these links will be 
identified as causal. We wish to investigate this smaller subset of links. 

Intuitively, if a time-series X causally affects another Y, then the past values of X 
should be helpful in predicting the future values of Y, above what can be predicted based 
on the past values of Y alone. This is formalized as follows: the error in predicting yt+s 
from yt, j/t-i, . . . should be larger than the error in predicting yt+s from the joint data 
yt,yti, ■ ■ ■ ,xt, Xt^i, ... if X causes Y. As our model, we chose to use the standard vector- 
autoregressive (VAR) model from econometrics [Box et al., 2011]. This is implemented 
in a Matlab code by LeSage [1999]. The standard p-lag VAR model takes the form: 



yt 

Xt 



Em. 



yt-i 

Xt-i 



where c is a vector of constants, Mi are the n x n coefficient (or autoregressive) mixing 
matrices and et is the unobservable white-noise. For the results shown below, p = 2. We 
then use the standard F-test to determine significance. 

In Table 4, we show the causal relationships identified among the out-links of the 
article Earthquake. Recall that there was a major earthquake in Australia during our 
time- window. We wish to understand which of the out-links appeared to be sensitive to 
this large change in interest in Earthquake. We use a significance cutoff of 0.01 and test 
for Granger causality among the time-series with s = 4. 



7 Conclusion 

PageRank is one of the most widely used network centrality measures. Our dynamical 
system reformulation of PageRank permits us to incorporate time-dependent teleporta- 
tion in a relatively seamless manner. Based on the results presented here, we believe 
this is an interesting variation on the PageRank model. For instance, we can analyze 
certain choices of oscillating teleportation functions (Lemma 6). Our empirical results 
show that the maximum change in the transient rank values identifies interesting sets of 
pages. Furthermore, this method is simple to implement in an online setting using either 
a forward Euler or Runge-Kutta integrator for the dynamical system. We hope that it 
might find a use in online monitoring systems. 

One important direction for future work is to treat the inverse problem. That is, 
suppose that the observed page views reflect the behavior of these random surfers. For- 
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Table 4 - Example of Causality in Wikipedia. We only consider pages with a pval < 0.01 as 
statistically significant. The page with values "caused" by Earthquake represent ideas related 
to earthquakes. All pages below are significant with pval < 0.01. 





n—T/Q 1 1 1 o 
Vdl Lit; 


Seismic hazard 


0.003535 


Extensional tectonics 


0.003033 


Landslide dam 


0.002406 


Earthquake preparedness 


0.001157 


Richter magnitude scale 


0.000584 


Fault (geology) 


0.000437 


Aseismic creep 


0.000419 


Seismometer 


0.000284 


Epicenter 


0.000020 


Seismology 


0.000001 



mally, suppose that we equate page views with samples of x(i). Then, the goal would 
be to find v(t) that produces this x(i). This may not be a problem for websites such 
as Wikipedia, due to our argument that the majority of page views reflect search engine 
traffic. But for many other cases, we suspect that x{t) may be much easier to observe. 
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